import schimpy
from schimpy import station
dfs = station.read_station_in('../tests/data/m1_hello_schism/station.in')
dfs = dfs.reset_index()
from schimpy import schism_mesh
grid = schism_mesh.read_mesh('../tests/data/m1_hello_schism/hgrid.gr3')
import holoviews as hv
hv.extension('bokeh')
from holoviews import opts, dim, streams
from holoviews.operation import datashader
import warnings
warnings.filterwarnings('ignore')
import panel as pn
pn.extension()
trimesh = hv.TriMesh((grid.elems, grid.nodes))
# rasterize to view faster, zoom in to clarify features
img = datashader.rasterize(trimesh.edgepaths).opts(opts.Image(cmap=['darkblue'])).opts(width=800, height=400)
# spread image pixels to see mesh in a more bold style
elems_only = datashader.spread(img)
nodes_only = datashader.dynspread(datashader.rasterize(trimesh.nodes).opts(opts.Image(cmap=['blue'])), shape='circle', max_px=6)
full_mesh = elems_only*nodes_only
elems_only.opts(alpha=0.2)*hv.Points(dfs, kdims=['x','y'],
vdims=['z','id','subloc','name']).opts(color='red',
size=10,
tools=['hover'])
trimesh = hv.TriMesh((grid.elems, grid.nodes))
# rasterize to view faster, zoom in to clarify features
img = datashader.rasterize(trimesh.edgepaths).opts(opts.Image(cmap=['darkblue'])).opts(width=800, height=400)
# spread image pixels to see mesh in a more bold style
elems_only = datashader.spread(img)
nodes_only = datashader.dynspread(datashader.rasterize(trimesh.nodes).opts(opts.Image(cmap=['blue'])), shape='circle', max_px=6)
full_mesh = elems_only*nodes_only
elems_only.opts(alpha=0.2)*hv.Points(dfs, kdims=['x','y'],
vdims=['z','id','subloc','name']).opts(color='red',
size=10,
tools=['hover'])
# param.nml has the time in the OPT section, but no parser for .nml files found
from datetime import datetime
import pandas as pd
import hvplot.pandas
def read_and_plot(file, station_file, reftime):
df1 = station.read_staout(file, station_file, reftime)
df1.index.name='Time' # workaround for hvplot bug
return df1.hvplot()
reftime = datetime(2000,1,1)
plots = []
for index in range(1,10):
fpath = f'../tests/data/m1_hello_schism/outputs/staout_{index}'
station_file = '../tests/data/m1_hello_schism/station.in'
vartype = schimpy.station.station_variables[index-1]
try:
plot = read_and_plot(fpath, station_file, reftime).opts(ylabel=f'{vartype}')
plots.append(plot)
except:
pass
#print(f'No data for index: {index}: {vartype}')
hv.Layout(plots).cols(1).opts(shared_axes=False)
selectable_pts = hv.Points(dfs, kdims=['x','y'],
vdims=['z','id','subloc','name']).opts(color='red',
size=10,
tools=['hover','tap'])
select_station = streams.Selection1D(source=selectable_pts)
def read_station_data(station_id, vartype, file_path_prefix, reftime):
station_index = schimpy.station.station_variables.index(vartype)
print(station_index)
file = f'{file_path_prefix}{station_index+1}'
df = station.read_staout(file, station_file, reftime)
df.index.name='Time' # workaround for hvplot bug
selected_cols = df.columns[df.columns.str.contains(station_id)]
return df[list(selected_cols)]
def show_station_data(station_id, vartype, file_path_prefix, reftime):
df = read_station_data(station_id, vartype, file_path_prefix, reftime)
return df.hvplot()
#show_station_data('ocean', 'elev', '../tests/data/m1_hello_schism/outputs/staout_', reftime)
def show_station_data_for_index(index, vartype):
if not index:
return hv.Div('Select a point from the map')
sdfs = dfs.iloc[index]
ids = sdfs['id'].unique()
return hv.Overlay([show_station_data(id,
vartype,
'../tests/data/m1_hello_schism/outputs/staout_',
reftime) for id in ids])
vartype = pn.widgets.Select(options=station.station_variables)
pn.Column(elems_only.opts(alpha=0.2)*selectable_pts,
vartype,
pn.Row(pn.bind(show_station_data_for_index,
index = select_station.param.index,
vartype = vartype)))
import numpy as np
points = hv.Points([])
taps = hv.streams.Tap(source=points, x=np.nan, y=np.nan)
def location(x, y):
if np.isnan(x) or np.isnan(y):
element_id = np.nan
else:
element_id = grid.find_closest_elems([x,y])
return pn.pane.Str(f'Click at {x:.2f}, {y:.2f}\nclosest element {element_id}')
pn.Column(points*elems_only, pn.bind(location, x=taps.param.x, y= taps.param.y))